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Summary We performed a genome-wide association study for Warner-Bratzler shear force (WBSF), 

a measure of meat tenderness, by genotyping 3360 animals from five breeds with 
54 790 BovineSNP50 and 96 putative single-nucleotide polymorphisms (SNPs) within 
ji-calpain [HUGO nomenclature calpain 1, (mull) large subunit; CAPN1] and calpastatin 
(CAST). Within- and across-breed analyses estimated SNP allele substitution effects (ASEs) 
by genomic best linear unbiased prediction (GBLUP) and variance components by 
restricted maximum likelihood under an animal model incorporating a genomic relation- 
ship matrix. GBLUP estimates of ASEs from the across-breed analysis were moderately 
correlated (0.31-0.66) with those from the individual within-breed analyses, indicating 
that prediction equations for molecular estimates of breeding value developed from 
across-breed analyses should be effective for genomic selection within breeds. We identi- 
fied 79 genomic regions associated with WBSF in at least three breeds, but only eight 
were detected in all five breeds, suggesting that the within-breed analyses were under- 
powered, that different quantitative trait loci (QTL) underlie variation between breeds or 
that the BovineSNP50 SNP density is insufficient to detect common QTL among breeds. 
In the across-breed analysis, CAPN1 was followed by CAST as the most strongly associ- 
ated WBSF QTL genome-wide, and associations with both were detected in all five breeds. 
We show that none of the four commercialized CAST and CAPN1 SNP diagnostics are 
causal for associations with WBSF, and we putatively fine-map the CAPN1 causal muta- 
tion to a 4581-bp region. We estimate that variation in CAST and CAPN1 explains 1.02 
and 1.85% of the phenotypic variation in WBSF respectively. 

Keywords beef, Bos taurus taurus, calpain 1, (mu/I) large subunit, calpastatin, genome- 
wide association, haplotype, meat tenderness, quantitative trait loci, single-nucleotide 
polymorphisms, Warner-Bratzler shear force. 



Introduction 

Consumer assessment of beef quality, palatability and 
overall eating satisfaction is significantly influenced by 
tenderness (Huffman et al. 1996; Weston et al. 2002; 
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Moser et al. 2004; Smith et al. 2006), and consumers 
have indicated a willingness to pay a premium for 'guar- 
anteed tender' steak (Boleman et al. 1997; Mintert et al. 
2000; Miller et al. 2001; Platter et al. 2005). Inadequate 
tenderness has consistently been identified in National 
Beef Quality Audits as a priority quality challenge 
(Lorenzen et al. 1993; Roeber et al. 2000; Shook et al. 
2008) because consumers consider tenderness to be the 
single most important component of meat quality and will 
substitute protein sources motivated by their dissatisfac- 
tion from the purchase of a tough cut (Miller et al. 1995; 
McKenna et al. 2002). 

To address these concerns, researchers have identified 
quantitative trait loci (QTL) for Warner-Bratzler shear 
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force (WBSF) measurements on the longissimus dorsi mus- 
cle on chromosomes 2, 4, 5, 7, 10, 11, 15, 20, 25 and 
29 (Casas et al. 1998, 2000, 2001, 2003; Keele et al. 
1999; Rexroad et al. 2001; Alexander et al. 2007; Davis 
et al. 2008; Gutierrez-Gil et a/. 2008; Gill et a/. 2009, 
2010). However, from these reported QTL, DNA marker 
tests have been developed and commercialized only from 
calpastatin (CAST) on chromosome 7 and calpain 1, (mull) 
large subunit (CAPNT) on chromosome 29 (Page et al. 
2002, 2004; White et al. 2005; Casas et al. 2006; Van 
Eenennaam et al. 2007). While these commercialized mar- 
ker tests are predictive of tenderness in both Bos taurus 
taurus and B. t. indicus breeds, it appears that they are 
not causal for the detected associations with tenderness 
(Casas et al. 2003). However, the estimated genotypic 
associations estimated for these markers are large, with 
an average difference of 0.15 kg in WBSF between alter- 
nate homozygotes in independent studies involving several 
breeds (Casas et al. 2006; Morris et al. 2006; Van 
Eenennaam et al. 2007; Johnston & Graser 2010). While 
positional candidate genes on other chromosomes have 
been investigated (Rexroad et al. 2001; Stone et al. 2005), 
none have resulted in commercial tests. 

To assist beef breeders to make efficient and large 
changes in tenderness, DNA assays must be developed 
that can reliably predict the genetic variation in tender- 
ness without regard to the breed composition of an 
animal. To address this need, we genotyped 3360 animals 
representing 114 half-sib families produced by the 
American Angus Association (AAA), American Hereford 
Association (AHA), American Simmental Association 
(ASA), American International Charolais Association 
(AICA) and the North American Limousin Foundation 
(NALF) as part of the National Cattlemen's Beef Associa- 
tion (NCBA) sponsored Carcass Merit Project (CMP) to 
develop prediction equations for the implementation of 
genomic selection (Meuwissen et al. 2001) and to identify 
genomic regions associated with tenderness. This study 
reports genomic regions detected as being concordant 
across breeds, which putatively harbour candidate genes 
that influence tenderness and which could be targeted for 
the development of diagnostic assays. We also dissect 
variation within CAST and CAPN1 in order to identify the 
genomic regions most likely to harbour the causal 
variants influencing beef tenderness. 

Materials and methods 

Animals and phenotype 

A total of 3360 animals representing five of the breed 
associations participating in the NCBA-sponsored CMP 
were selected for genotyping based on the availability of 
WBSF data and DNA samples (Table 1). The design of the 
CMP project has previously been described by Minick et al. 
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(2004); however, only the Angus and Hereford samples 
represent purebred populations, with the Continental 
breeds being represented by crossbred progeny, with Sim- 
mental, Charolais and Limousin sires mated to predomi- 
nantly commercial Angus cows. Meat tenderness was 
measured as WBSF (kg) of longissimus dorsi steaks at day 
14 post-mortem as previously described (Wheeler et al. 
1998; Minick et al. 2004). Muscle samples, extracted 
DNA samples and carcass phenotypes produced in the 
CMP and owned by the AAA, AHA, ASA, AICA and 
NALF were transferred to the University of Missouri. All 
CMP animals had blood samples drawn at weaning, from 
which DNA was extracted and tested to validate the iden- 
tity of their sires. Additionally, a muscle sample was taken 
at slaughter at the capture of phenotype data on most of 
the animals, and DNA extracted from a subset of the mus- 
cle samples was previously genotyped and compared with 
the genotype profiles produced from the corresponding 
blood samples to validate the identity of each carcass. This 
process identified that about 10% of animals or carcasses 
were misidentified (Thallman et al. 2003) likely due to 
changes in the order of carcasses because of 'rail-outs' at 
packing plants. To resolve this issue, we extracted geno- 
mic DNA from 2940 muscle samples taken from the phe- 
notyped carcasses by proteinase K digestion followed by 
phenol-chloroform-isoamyl alcohol extraction and etha- 
nol precipitation (Sambrook et al. 1989). The remaining 
420 DNA samples were extracted from the blood, but 
these samples had previously been DNA-typed and 
successfully matched to the sample taken at harvest. 

Genotypes 

All samples were genotyped using the Illumina 
BovineSNP50 BeadArray (Matukumalli et al. 2009) for 
54 790 single-nucleotide polymorphisms (SNPs) and a cus- 
tom-designed Illumina GoldenGate assay incorporating 96 
putative SNPs located within 186 kb of CAST and CAPN1 
(Table SI), discovered either as part of the bovine genome 
sequencing project or through directed CAPN1 resequenc- 
ing studies at the US Meat Animal Research Center at Clay 



Table 1 Animal counts, mean phenotype and estimates of additive 
genetic variance and heritability by breed. 



Breed 


Count 




Warner-Bratzler shear 
force (kg) 




Animals 1 


Sires 


Average 


°l 


h 2 


Angus 


660 (651) 


20 


3.74 


0.22 


0.52 


Charolais 


702 (695) 


18 


4.41 


0.23 


0.46 


Hereford 


1192 (1095) 


29 


4.75 


0.15 


0.17 


Limousin 


285 (283) 


23 


4.28 


0.07 


0.09 


Simmental 


521 (516) 


24 


4.36 


0.06 


0.08 


All Breeds 


3360 (3240) 


114 


4.37 


0.17 


0.25 



1 Numbers of animals with genotype call rate >0.85 in parentheses. 
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Center, NE (Page et al. 2002; White et al. 2005; Casas et al. 
2006). Several of the putative SNPs identified in the gen- 
ome sequencing project were not variable (Table SI), and 
we were much more successful in fine-mapping CAPN1 
than CAST. All genotypes were called in the lllumina 
genomestudio software. Genotypes were filtered according to 
their unique localization to an autosome or the X chromo- 
some in the University of Maryland sequence assembly 
(UMD3.0; Zimin et al. 2009), call rate (>0.89) and minor 
allele frequency >0.01 within each breed. Animals were 
excluded if their individual genotype call rate was <0.85. 
The call rate of >0.89 for SNP filtering was used to ensure 
that all commercialized tenderness SNPs were included in 
the analysis. After filtering, the data set comprised 40 645 
SNPs assayed in 3240 animals (Tables 1 and S2). 

Analysis 

fastphase vl.2.3 (Scheet & Stephens 2006) was used with 
UMD3.0 coordinates to phase all genotypes and impute the 
0.89% of missing genotypes. The complete set of genotypes 
was then used to generate a genomic relationship matrix 
(G) across all breeds using the first of the methods 
described by VanRaden (2008) with a modification allow- 
ing the inclusion of X-linked loci as described below. 

Warner-Bratzler shear force phenotypes were analysed 
under a single-trait mixed linear animal model in which 
the genomic relationship matrix was used to represent the 
realized identity by descent among the animals. The model 
fit was y = Xfi + Zu + e where y is a vector of WBSF mea- 
surements, \i is a vector of fixed contemporary group effects 
defined as breed x herd of origin x sex of calf x slaughter 
date, u is a vector of random additive genetic merits, 
and e is a vector of random residuals. The matrices X 
and Z are incidence matrices relating observations to lev- 
els of the fixed and random effects, and we assume that 
Var(u) = Gff^,Var(e) = Iffg andCov(u,e) = 0. Restricted 
maximum likelihood was used to estimate the variance 
components a\ and a\ and iteration on the variance com- 
ponent estimates continued until the estimate of heritabil- 
ity h 2 = a\/(a\ + a 2 ,) had converged to four significant 
figures. At convergence, the GBLUP of the vector of SNP 
allele substitution effects (ASEs) was obtained as 
a. = (2SiP i q i ) _ M'G~ ii where p t is the frequency of the A 
allele at the ith SNP (genotypes at each SNP are called in 
A/B space by the GenomeStudio software), <j f = 1 - p b ele- 
ments of the ith column of M are 2q h q t — p t and — 2p t for 
AA. AB and BB genotypes at autosomal and pseudoautos- 
omal loci (VanRaden 2008) and are q t and -p t for AY and 
BY genotypes at X-linked loci in males, and u is GBLUP of 
u. Analyses were performed both within each breed and 
across all breeds. 

The variance component associated with SNP ASEs is 
aju = (2Z i p i q i ) _1 f7^, and for each SNP, the predicted ASE 
was normalized to a t-like statistic as t t = | at; | la M . These 



values are included in Table S2 and are shown in the 
Manhattan plots in Figs 1 and SI. 

Across-breed comparison of putative QTL regions 

To determine whether common QTL influence WBSF across 
breeds, we ranked the t t values estimated in the within- and 
across-breed analyses and then identified SNPs for which 
the t,- values ranked in the top 500 (1.2%) of SNP ASEs in 
the across-breed analysis. For each of the regions tagged by 
these SNPs, we declared the region to harbour a QTL if at 
least three SNPs from different within-breed analyses had 
ASEs ranked in the top 500. While linkage disequilibrium 
(LD) decays to -0.1 within less than a 500-kb distance 
within breeds of distantly related individuals (McKay et al. 
2007), many of the individuals incorporated into these 
analyses are half-sibs (Table 1), which leads to a much 
greater extent of LD because of large common chromosomal 
segments transmitted by the sires to their progeny. Addi- 
tionally, we wanted to allow for the possibility that more 
than one QTL could be present within any one genomic 
region. Accordingly, we allowed the region size to vary up 
to 5.7 Mb (average 1.7 Mb) as determined by the signa- 
tures of the detected within-breed SNP ASE ranks. Further- 
more, within each region, we did not expect to find the 
same SNP to be most strongly associated with WBSF, 
because differences in SNP and QTL allele frequencies 
between breeds (Table S2) can lead to different patterns of 
LD in different breeds. 

Candidate genes 

Genomic regions identified as being associated with WBSF 
in at least four breeds were analysed using the NCBI 

1.00 n 
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0.80 - 
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Cumulative genome position (Mb) 

Figure 1 Manhattan plot of single-nucleotide polymorphism 
(SNP) allele substitution effects estimated in the across-breed 
analysis and normalized by the square root of the estimated SNP 
variance component. 
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Entrez Map Viewer (accessed 07/06/2011) to identify 
potential candidate genes for tenderness. 

CAST and CAPN1 

A 1.48-Mb region of BTA7 harbouring 28 SNPs spanning 
CAST and a 2.64-Mb region of BTA29 harbouring 93 
SNPs spanning CAPN1 were found to contain loci for 
which SNP ASEs ranked in the top 500 in the within- 
breed analyses. To allow haplotype-based analyses, we 
expanded the regions to 44 SNPs spanning 2.86 Mb for 
CAST and 100 SNPs spanning 3.12 Mb for CAPN1 (Table 
S3). We first analysed each SNP individually by including 
allele effects (the difference between the two estimated 
allele effects is the ASE for the SNP) in ft, in addition to 
the contemporary group effects, and then we included 
haplotype effects for windows of nine contiguous SNPs 
using phase information estimated by fastphase. The hap- 
lotype model was sequentially fit by sliding the nine SNP 
window through each region one SNP at a time, and the 
statistics computed for each window were assigned to the 
5th SNP located at the centre of each window. In both 
cases, the analysis was performed using the previously 
estimated variance components (Table 1), and F-tests for 
SNP or haplotype effects were constructed from the differ- 
ence between model sums of squares including and 
excluding the fitted SNP or haplotype effects, the difference 
in number of parameters between the fitted models and 
the estimated residual variance for the full model. Because 
the number of detected haplotypes varied throughout each 
region (Table S3), the window producing the largest 
model sum of squares does not necessarily result in the 
largest F-statistic or — log 10 P-value (because the numera- 
tor mean square can be significantly influenced when its 
degrees of freedom are small but vary between tests). To 
avoid this, we computed the percentage of phenotypic var- 
iation explained by each window through the region from 
the ratio of the window to phenotypic sums of squares, 
where the window sum of squares was estimated as the 
difference between model sum of squares including and 
excluding haplotype effects for the nine SNP window and 
the phenotypic sum of squares was estimated as the total 
sum of squares corrected for the mean and contemporary 
group sums of squares. This statistic identifies the SNP 
window that explains the largest amount of variation in 
WBSF regardless of the number of haplotypes that are fit. 

Results and discussion 

We found large differences in the heritabilities of WBSF 
across the five breeds (Table 1) and were concerned that 
this might reflect differences in data quality or the correct 
assignment of phenotype to genotype because of the sam- 
ple misidentification issue identified within the CMP. How- 
ever, we also estimated heritabilities for eight additional 



carcass traits recorded in this project (data not shown) 
and found no evidence for systematically lower heritabili- 
ties within any of the breeds. We therefore conclude that 
the re-extraction of DNA from tissue samples taken from 
the carcass at slaughter effectively solved the misidentifi- 
cation problem. Thus, the variation in heritabilities proba- 
bly reflects the relatively small sample size within each 
breed and the sampling of the bulls used to produce these 
animals. However, the effect of variation in heritability 
across breeds was to substantially influence the 'genetic' 
sample size which we estimate as N x h, the number of 
phenotypes multiplied by the square root of the heritabil- 
ity, which is an estimate of the cumulative amount of 
additive genetic information in a sample of N unrelated 
individuals and was 468.3, 451.5, 471.7, 85.4 and 
143.5 in Angus, Hereford, Charolais, Limousin and Sim- 
mental respectively. 

In the across-breed analysis, the use of the genomic 
relationship matrix corrects for the stratification because 
of pedigree relatedness while accounting for the extent of 
background relatedness among the Angus and Continen- 
tal breed groups because of the use of Angus dams to pro- 
duce the crossbred Continental breed calves. In this 
analysis, the associations between the CAST and CAPN1 
loci with WBSF were the largest in the genome (Fig. 1), 
reflecting both the magnitude of effects of these genes and 
the increased SNP density within these regions, which 
improves the likelihood of finding SNP in strong LD with 
the causal mutations. The within-breed analyses identified 
CAPN1 as the locus most strongly associated with WBSF 
genome-wide, although the highest ranked SNP ASE 
within this region for Limousin was only 30th (Table S2), 
presumably reflecting the very small sample size for this 
breed. On the other hand, the CAST associations were 
more variable among the breeds, being the most strongly 
associated with WBSF genome-wide in Hereford, ranking 
highly in Charolais and Limousin, but only 234th and 
208th in Angus and Simmental respectively. These results 
are likely due to the fairly small sample sizes for the anal- 
ysed breeds, but probably also may reflect the different 
SNP densities within the two regions and differences in 
allele frequencies at the SNPs and QTL across breeds. We 
accomplished a much higher SNP density in the region 
harbouring CAPN1 than CAST, and this suggests that we 
had insufficient SNPs to find at least one that was in 
strong LD with the causal mutations within CAST in all 
breeds. 

Across all 40 645 SNPs, the correlations between ASEs 
estimated within each of the breeds varied from —0.02 to 
0.04, indicating that models developed to predict genomic 
breeding values within one breed will have very low accu- 
racies in other breeds. This has previously been predicted 
using simulated data (de Roos et al. 2009; Toosi et al. 
2010) but, despite the use of commercial Angus females 
to produce the Continental breed crossbred steers, it is a 



© 2012 The Authors, Animal Genetics © 2012 Stichting International Foundation for Animal Genetics, 43, 662-673 



666 McClure et al. 



consequence of the genetic distance between the training 
and validation sets of animals. Habier et al. (2010) dem- 
onstrated that the number of generations that separate 
the training and validation data sets influences the 
accuracy of genomic breeding values estimated in the 
validation set, with lower accuracies occurring when this 
relationship is more distant. On the other hand, the corre- 
lations between the ASEs estimated in the across-breed 
analysis and those estimated in the within-breed analyses 
were 0.37, 0.66, 0.41, 0.31 and 0.42 for Angus, 
Hereford, Charolais, Limousin and Simmental respectively. 
This result supports the simulation results of Toosi et al. 
(2010), who showed that training in admixed populations 
results in genomic estimates of breeding value with 
accuracies almost equivalent to those achieved from 
training and validating within the same breed. Of course, 
the key benefits from the perspective of beef cattle breed- 
ing are that training population samples can dramatically 
be increased by pooling breeds and that the resulting 
genomic breeding values have industry-wide utility. 

Hayes & Goddard (2001) have estimated that between 
50 and 100 QTL underlie variation in quantitative traits 
within livestock populations. While under neutral theory, 
the common QTL mutations that are detectable by GWA 
analysis must predate the domestication of cattle (Kimura 
& Ohta 1973), the relatively small populations upon 
which breeds were founded may have led to the sampling 
of different subsets of QTL within different breeds. In fact, 
the extent to which breeds share common QTL is 
unknown (Pryce et al. 2010), but is of some importance 
to the development of prediction equations for molecular 
estimates of breeding value in admixed populations and 
the development and utilization of genotyping assays for 
the prediction of genetic merit within the beef industry. 
To identify QTL underlying variation in WBSF, we exam- 
ined the genomic regions harbouring the 500 SNPs with 
the largest ASEs from the across-breed analysis for SNPs 
with ASEs ranked in the top 500 in the within-breed anal- 
yses for at least three of the breeds. We identified 79 geno- 
mic regions that putatively harbour QTL influencing 
WBSF (Table 2). Of these, 42 were identified in three 
breeds, 29 in four breeds and eight in all five breeds. 
There was no difference between the breeds (P = 0.48) or 
between British and Continental breeds (P = 0.52) in the 
probability of QTL detection for all 79 QTL or for the 42 
QTL identified in only three breeds (P = 0.35 and 0.82 
respectively). Clearly sample size, assay SNP density, con- 
straints on SNP ranks and the size of regions harbouring 
highly ranked SNP ASEs all impact the identification of 
putatively common QTL. Of the 113 instances when the 
within-breed estimated SNP ASEs ranked >500, the aver- 
age rank was only 2551, suggesting that the majority of 
these regions harbour QTL that segregate in all breeds. 
Changing the minimum within-breed ASE rank criterion 
to <1000 resulted in 17 of these QTL being detected in all 



five breeds, 41 in four breeds and 21 in three breeds 
(Table 2). Thus, there appears to be little phylogenetic sig- 
nal in these data, and if a QTL was detected in only three 
breeds, these breeds were as likely to be British and Conti- 
nental as strictly Continental. 

We have previously found poor concordance between 
GWA and half-sib linkage analyses for large-effect QTL 
underlying growth traits, even when large numbers (>50) 
of families with family sizes ranging from 20 to 224 half- 
sibs are analysed (data not shown). Assuming that GWA 
analysis detects common variants, we would expect a sig- 
nificant number of sires to be both heterozygous and 
detected to be segregating for a large-effect QTL; however, 
this largely depends on the underlying genetic architec- 
ture of the trait. Reed et al. (2008) found that growth was 
affected in 34% of viable mouse knockouts, suggesting 
that natural variation in thousands of genes underlies var- 
iation in growth. As a consequence of this complex 
genetic architecture, there may be a large number of QTL 
on each chromosome, and the allelic combinations present 
at these QTL in the sire will impact on whether any one 
QTL is detected in linkage analyses. Thus, common vari- 
ants detected in GWA analysis may not be detected in seg- 
regation analysis, and rare variants detected in 
segregation analysis may not be detected in GWA analy- 
sis. Nevertheless, we found six of the 12 previously 
reported meat tenderness QTL, including CAST and 
CAPN1, to coincide with the QTL identified in this 
study (Table 2) (Cattle QTL database, http://www.animal 
genome.org/ cgi-bin/ OTLdb/BT/ draw_traitmap?trait_iD= 
1030, accessed June 27, 2011). Notwithstanding the poor 
resolution of QTL location mapped by linkage analysis, we 
also found support for all of the other previously identified 
QTL. For example, in the across-breed analysis, QTL were 
identified with ASE ranks < 5 00 at 3 151 989 bp and at 
6 831 955-7 086 105 bp (300 kb from MSTN) on 
BTA2. The first was supported by ASE ranks <500 for 
Angus and Charolais, but an ASE rank of 565 in 
Limousin. The second was supported by an ASE rank 
<500 in Charolais and ASE ranks <1000 in Angus, Lim- 
ousin and Simmental. Thus, despite their proximity, these 
QTL are likely distinct, and the concordance between our 
and previously published results suggests that the genetic 
architecture of meat tenderness is substantially less 
complex than for growth. 

We examined the genomic regions harbouring the 37 
QTL that were detected in at least four of the breeds for 
potential candidate genes underlying meat tenderness. 
Very little is known about the genetic regulation of meat 
tenderness, and few candidate genes are suggested for 
these QTL. While CAST and CAPN2 have consistently 
been identified and analysed as candidate genes for the 
BTA7 97 861 341-98 820 742-bp and BTA29 
44 042 363-44 087 629-bp QTL, respectively, no causal 
variants have been identified in either gene. CAPN1 
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encodes the protease /i-calpain, which has been implicated 
in the proteolysis of muscle proteins during meat ageing 
(Smith et al. 2000), and CAST encodes calpastatin, which 
is an inhibitor of ^i-calpain (Goll et al. 2003). Myogenic 
determination factor 1 is a transcription factor encoded by 
MY0D1 and is expressed in skeletal muscle during myo- 
genesis and regeneration. Variation in MY0D1 has been 
suggested to affect its ability to influence the expression of 
muscle structural components (Rexroad et ah 2001), mak- 
ing it a candidate for the QTL at 34682 617- 
36 817 688 bp on BTA15. Calpain-2 (mill) large subunit 
(m-calpain) is a calcium-activated neutral protease 
encoded by CAPN2 on BTA16 (25 000 153- 
28 384 914 bp). M-calpain activity has been associated 
with both meat tenderness and palatability measurements 
(Riley et al. 2003). Fibroblast growth factor 2 (FGF2) is 
an upstream regulator of heat shock protein Bl (HSPB1), 
which has been found to be negatively related to WBSF 
(Kim et al. 2011), making it a candidate for the 
34 429 947-37 201 424-bp QTL on BTA17. GSN 
encodes gelsolin, a calcium-regulated protein that func- 
tions in both the assembly and disassembly of actin fila- 
ments, which are a component of the contractile 
apparatus in muscle cells and may underlie the BTA8 
112 287 843-113 301 368-bp QTL. Finally, CALM1 
encodes calmodulin, a calcium-binding protein, which 
interacts with titin and mediates smooth muscle contrac- 
tion, making it a candidate for the BTA10 102 286 251- 
103 234 411-bp QTL. 

While the commercially tested CAST SNP rs41255587 
was the most strongly associated with WBSF in the 
across-breed analysis (— log 10 P = 8.95), it was only the 
most strongly associated CAST SNP within Hereford and 
Charolais. with stronger associations being detected for 
SNPs in the 5' upstream region in Angus. Limousin and 
Simmental (Table S3). In fact, the haplotype analysis 
moves the location of the most significantly associated 
SNP window 83.7 kb upstream of rs4 125 5 587 to be cen- 
tred on rs43529872 (-log 10 P = 8.78), and this CAST 
window was found to explain the greatest amount of phe- 
notypic variation in WBSF in the across-breed (1.02%; 
Table 3), Angus and Hereford analyses. The sign and 
magnitude of the ASE was consistent for rs41255587 in 
all breeds except Limousin, and the haplotype analysis 
explained considerably more variation in WBSF than the 
single SNP analysis, indicating that either the causal vari- 
ant is not among the tested polymorphisms or that there 
is more than one causal variant. Furthermore, the haplo- 
type analyses move the most likely location of the causal 
mutation 5' of the commercially tested CAST SNP 
rs41255587, probably in the 678-kb region from 
97 861 341-98 538 952 bp (Fig. 2). Clearly, additional 
fine-mapping is required to identify the number of muta- 
tions influencing WBSF that lie in the vicinity of CAST 
and their most likely locations. 
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Figure 2 Proportion of phenotypic variation in the across-breed 
analysis explained by haplotypes constructed from nine consecutive 
single-nucleotide polymorphism (SNPs) in the region of (a) BTA7 
harbouring CAST and (b) BTA29 harbouring CAPN1. tocations and 
amount of variation explained by the commercialized tenderness 
SNPs are indicated by red dotted lines. 

Among the SNP located within CAPN1. rsl7812000 
(c.316G>A) was most strongly associated with WBSF in 
Angus (-logioP = 9.70) and rsl7872050 was the most 
strongly associated with WBSF in Limousin (— log 1() P = 
3.23). However, rs42192103 was found to be slightly 
more strongly associated with WBSF than rsl7 812000 in 
the across-breed analysis (— log 10 P = 15.25 vs. 15.01), 
with an average ASE across breeds of 0.2 3 kg (Table S3). 
The amount of phenotypic variation explained in 
the haplotype-based analyses again indicates that none of 
the tested SNPs are causal for effects on WBSF and that 
the strongest signal for association with WBSF was in the 
8187-bp region from 44 062 694 to 44 070 881 in all 
five breeds (Table S3). The size of this region is sufficiently 
small to speculate that there is probably only a single 
mutation in CAPN1 affecting WBSF in all Bos t. taunts 
cattle breeds, and the across-breed haplotype analysis 
shown in Table S3 and Fig. 2 suggests that the most 
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Table 3 Percentages of phenotypic variation in WBSF explained by the commercialized SNPs, the most strongly associated SNPs and haplotypes 
within the most strongly associated nine SNP window within CAST and CAPN1. 



Locus 


All I I ,. 

All breeds 


Angus 


Hereford 


Charolais 


Limousin 


Simmental 


CAST (BTA7) 














rs472555S7 1 


0.66 


0.53 


1 .47 


1.14 


0.70 


0.02 


98 579 574 














SNP 2 


0.66 


0.54 


1.47 


1.14 


2.28 


1.13 




98 579 574 


98 498 047 


98 579 574 


98 579 574 


97 861 341 


98 013 150 


Window-P 3 


1.02 


1.36 


1.88 


2.10 


3.88 


2.77 




98 495 888 


98 495 888 


98 566 391 


98 538 952 


97 501 859 


97 861 341 


Window-Vp 4 


1.02 


1.36 


1.92 


2.10 


4.02 


2.77 




98 495 888 


98 495 888 


98 495 888 


98 538 952 


98 375 640 


97 861 341 


CAPN1 (BTA29) 














rs17812000 1 


1.14 


2.36 


0.96 


1.38 


0.00 


3.75 
















rs1787105l' 1 


0.39 


1.54 


0.16 


0.39 


0.57 


1.66 


44 085 642 














rs1 7872050" 


0.53 


0.89 


0.08 


1.21 


2.88 


1.65 


44 097 629 














SNP 2 


1.16 


2.36 


1.62 


1.57 


2.88 


4.65 




44 070 713 


44 069 063 


44 067 796 


44 070 713 


44 087 629 


44 042 363 


Window-P 3 


1.80 


3.18 


2.59 


2.76 


2.99 


5.05 




44 067 796 


44 068 519 


44 062 694 


44 070 881 


44 087 356 


44 067 234 


Window-Vp 4 


1.85 


3.19 


2.59 


2.76 


3.52 


5.35 




44 068 143 


44 068 445 


44 062 694 


44 070 881 


44 070 881 


44 068 143 



CAST, calpastatin; CAPN1 , calpain 1, (mu/l) large subunit; SNP, single-nucleotide polymorphism; WBSF, Warner-Bratzler shear force. 

Commercialized SNP and its chromosomal coordinate. 

2 Most strongly associated SNP and its chromosomal coordinate. 

3 Most strongly associated nine SNP window centred on SNP with shown chromosomal coordinate. 
4 Nine SNP window explaining the greatest amount of phenotypic variation in WBSF. 



likely region harbouring this mutation is the 4581-bp 
region from 44 063 938 to 44 068 519. This region is 
wholly contained within the 5' end of CAPN1. We 
estimate from the haplotype analysis that CAPN1 explains 
1.85% of the phenotypic variation in WBSF in taurine 
cattle (Table 3). 

Conclusions 

We conclusively demonstrate that none of the SNPs cur- 
rently commercialized as diagnostics for genetic merit are 
causal for their effects on WBSF (Casas et al. 2003, 2006; 
Van Eenennaam et al. 2007; Gill et al. 2009). In fact, the 
complex patterns of LD in the vicinity of these genes 
among the different breeds (Figs S2 and S3) and the 
weaker associations in Limousin and Simmental (Fig. SI) 
result in different SNPs being most strongly associated 
with WBSF among the breeds (Table S3). However, by 
using haplotype-based analysis methods to dissect the var- 
iation within these genes, we localized the causal variants 
to be 5' to the commercially tested SNPs. In the case of 
CAPN1, the higher SNP density achieved and the use of 
across-breed analysis, which erodes the patterns of LD 
within breeds, resolved the likely location of the causal 
variant to a region of only 4581 bp. 



We found evidence for a large number of QTL underly- 
ing variation in WBSF, and the majority of the previously 
published QTL were validated in this analysis. We found 
reasonably strong evidence that most QTL were segregat- 
ing in all five breeds; however, the small genetic sample 
sizes for Limousin and Simmental make this comparison 
problematic, and it remains an unanswered question as to 
the extent to which breeds may share private alleles at 
QTL. This has previously been found in Belgian Blue, 
Marchigiana and Piedmontese cattle, where breed-specific 
polymorphisms in MSTN produce the double muscled phe- 
notype (Grobet et al. 1997; Kambadur et al. 1997; 
McPherron & Lee 1997; Marchitelli et al. 2003). This 
issue is of importance to the development of prediction 
equations for molecular breeding values in across-breed 
analyses, because the ASEs estimated for QTL regions will 
be averaged across breeds that segregate and those that 
do not segregate for certain QTL, which will limit the 
accuracy of molecular estimates of breeding value. Despite 
this, we found moderate correlations between GBLUP pre- 
dictions of ASEs computed in the across- and within-breed 
analyses, suggesting that the BovineSNP50 assay has suf- 
ficient resolution for the development of prediction equa- 
tions for genomic selection in beef cattle despite their 
considerably larger effective population size relative to 
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dairy cattle (The Bovine HapMap Consortium 2009), and 
also that WBSF QTL are commonly shared among breeds. 

Despite the apparent reduced complexity of a trait such 
as meat tenderness relative to growth, there appear to be 
a large number of QTL underlying variation in WBSF, and 
the identification of all of the mutations that underlie 
these QTL might appear to be intractable. However, recent 
developments in high-density SNP genotyping, high- 
throughput sequencing and genotype imputation suggest 
new strategies for the rapid simultaneous identification of 
variants underlying quantitative traits genome-wide. We 
accomplished an average SNP spacing of 1139 bp for the 
23 SNPs analysed within CAPN1, and this is only slightly 
smaller than could be accomplished genome-wide by 
jointly genotyping with the newly available Illumina Bovi- 
neHD and Affymetrix BOS 1 assays (-1.3 million SNP, 
data not shown). Furthermore, the design of these assays 
was facilitated by a community effort that produced more 
than 128. 4X of genome sequence coverage on more than 
80 animals, and SNP data from this work are now avail- 
able in dbSNP. This project discovered 48.6 million high- 
quality SNPs, which must include many of the causal 
variants underlying quantitative variation in cattle, and it 
may be possible to impute genotypes at the resolution of 
the genome sequence (Daetwyler et al. 2011) in popula- 
tions that have been genotyped with both assays. Such a 
strategy could rapidly allow the identification of a large 
number of causal variants if the association analysis was 
performed in mixed breed populations. 
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